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Computer  Simulations  of  Ion  Transport  across 
Excitable  Membrane  Stimulated  by  External  Field 


This  is  the  Final  Technical  Report  on  work  performed  under  the  support  of  the  Air 
Force  Office  of  Scientific  Research  (F49620-03-0134)  on  the  subject  that  was  changed  to 
Terahertz  radiation  sources  after  fulfillment  of  its  original  biological  direction.  This 
report  describes  the  motivation  and  approach  being  taken  to  elucidate  the  physical 
mechanism  for  generating  radiation  from  DC  discharge  in  Coaxial  Resonators.  The 
electrostatic  phenomena  arose  in  this  configuration,  such  as  plasma  oscillations  and 
sheath  formation  was  extensively  investigated  by  Irving  Langmuir1,2.  Hopefully,  we  will 
be  able  to  study  the  electromagnetic  effect  using  our  particle  model  under  development. 


I.  Introduction 

Terahertz  radiation  is  loosely  defined  in  the  frequency  range  of  0.3  to  3  THz 
(wavelength  range:  1000~100p,m).  Being  situated  between  infrared  light  and  microwave 
radiation,  THz  source  is  resistant  to  the  generation  and  diagnostic  techniques  commonly 
employed  in  these  well-established  neighboring  bands.  The  most  significant  limitation  to 
the  application  of  modem  THz  systems  is  the  lack  of  high-power,  low-cost,  room- 
temperature  compact  THz  source.  One  of  the  primary  motivations  for  the  development 
of  THz  source  is  its  potential  to  extract  material  characteristics  that  are  unavailable  when 
using  other  frequency  band.  In  recent  years,  THz  spectroscopy  systems  have  been 
applied  to  a  large  variety  of  materials  both  to  aid  the  basic  understanding  of  the  material 
properties,  and  to  demonstrate  potential  applications  in  sensing  and  diagnostics.  Other 
applications  include  Terahertz  tomographic  imaging  to  replace  mammography,  screening 
passengers  for  explosive  at  airports,  and  broadband  cell-phone  communication. 

Currently,  there  is  a  vast  array  of  potential  sources3  each  with  its  own  advantages. 
Broadband  THz  radiation  can  be  produced  by  thermal  sources  and,  more  recently,  by 
table-top  laser-driven  sources  using  high-speed  photoconductors  as  transient  current 
sources  for  radiating  antennas.  Using  these  approaches,  both  attainable  output  power 
(nano  to  microwatt  range)  and  efficiency  (10'6)  are  very  low.  In  the  last  decade,  active 
research  has  been  devoted  to  produce  narrow  band,  high  power  (kW),  and  high  efficiency 
(one  percent)  based  on  free-electron  laser  mechanism4.  In  this  device,  the  acceleration  of 
the  electron  beam  (2MV)  is  provided  by  a  small-size  RF  LINAC.  The  experimental  set¬ 
up  can  be  placed  on  a  bench-top  and  the  total  area  occupied  is  about  1x2  m2. 

An  alternative  approach  to  generate  THz  radiation  was  proposed5  in  1985  by 
Alexeff.  He  claimed  that  he  obtained  one  watt  output  with  10‘5  efficiency  Terahertz 
radiation  from  his  table-top  DC  (lOkV)  discharge  in  a  coaxial  resonator  without  magnetic 
field  (Fig.l).  However,  the  physical  mechanism  responsible  for  generating  the  THz 
radiation  in  this  configuration  still  remains  controversially.  There  are  the  following  three 
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as  it  transverses  the  sheath  region  around  the 
anode  wire7.  In  microwave  range,  the 
performance  based  on  DC  discharge  is 
unable  to  compete  against  the  conventional 
devices  in  output  power  and  efficiency. 
However,  most  experimental  results 
published  were  in  this  region.  In  order  to  be 
able  to  validate  simulation  results, 
simulations  will  be  carried  out  in  GHz  range 
first  to  isolate  the  physical  mechanism.  The 
parameters  will  then  be  scaled  up  to 
simulate  THz  generation. 


II.  Review  of  Experimental  Results 

Three  experimental  results  that  are  relevant  to  the  electromagnetic  radiation 
generation  based  on  the  concentric  cylindrical  diode  configuration  will  be  reviewed  here. 
A  DC  discharge  across  a  diode  can  be  pure  electron  discharge  or  discharge  in  gases. 


1 .  Orbitron 


The  basic  principle  of  operation  of  an  orbitron  is  believed  to  be  that  electrons 
oscillate  (orbital  motion)  in  a  system  in  which  the  electrons  both  are  confined  for  a  long 
time  and  the  oscillation  frequency  varies  with  energy.  When  an  electron  interaction  with 
the  cavity  mode  of  the  surrounding  metal  container  to  give  energy  to  the  cavity  EM  field, 
it  will  slow  down  and  then  move  closer  to  the  center  conductor  since  the  electrostatic 
field  force  will  be  greater  than  the  centrifugal  force.  The  electron  will  then  have  a 
slightly  smaller  equilibrium  radius  and  corresponding  higher  angular  velocity.  It  can  now 
further  interact  with  the  cavity  EM  field  in  approximate  synchronism  with  it  and  this 
process  may  continue  until  the  electron  is  intercepted  by  the  center  conductor.  By 
balancing  the  centrifugal  force  and  the  electrostatic  force  in  a  coaxial  system,  we  have. 


mv2 

r 


=  -eEr  =  - 


eV0 

rln(^) 

*0 


(1) 


where  ri  is  the  radius  of  the  cavity  and  r0  is  the  radius  of  the  wire.  The  electron  orbital 
frequency  is  given  by 
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From  Eq.(2),  the  highest  frequency  that  can  be 
obtained  is  when  r  =  ro.  The  original  orbitron 
design  is  shown  in  Figure  1.  By  using  the 
following  parameters:  ro=  8.5xl0'6m,  Vo=10kV, 
and  P(gas  pressure)  =  2x10  mm,  the  resulting 
current  is  20A  and  device  output  power  is  about 
1  watt  (eff.  =  10‘5)  at  1  THz. 
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Using  ri=  5x10’  m,  the  highest  frequency 

becomes  f0  =  6.5x1 0n  Hz.  On  the  other  hand,  Fig.  2.  Spectrum  of  gas  filled  steady 
the  plasma  frequency  for  a  fully  ionized  static  Orbitron. 

background  gas  is  about  2.7x1 011  Hz.  The 
fourth  harmonic  of  plasma  oscillation  can  also 

reach  1  THz.  In  a  later  experiment  ,  the  wire  radius  is  increased  to  2.5xl0‘  m  and  the 
applied  voltage  is  reduced  to  500V.  In  doing  so,  the  current  drops  to  only  a  few  mA  and 
the  orbitron  emits  a  surprisingly  stable  line  spectrum  of  a  few  GHz,  as  shown  in  Fig.  2. 

The  observed  lowest  frequency  mode  is  believed  to  be  due  to  coupling  to  the  external 
circuit.  The  discharge  current  is  observed  to  affect  the  output  frequency.  However, 
experimental  results  do  not  provide  any  conclusive  evidence  that  the  physical  mechanism 


for  radiation  generation  is  due  to  electron  orbital 
motion  around  the  positively  charged  wire. 

2.  Diverging  Concentric  Diodes 

The  theoretical  effort  to  produce  negative 
impedance  by  means  of  an  electron  discharge  in  a 
plane  diode9  has  been  known  since  1928.  The 
negative  resistance  appears  whenever  the  electron 
transit  time  between  electrodes  is  approximately 

I4 , 2^  ,  3^  etc  cycles  of  a  given  high-frequency 

current.  This  concept  was  applied  to  a  cylindrical 
configuration10  (Fig.  3  a)  in  which  the  inner 
conductor  is  the  emitter  (cathode)  while  the  outer 
conductor  is  the  collector  (anode)  that  is  biased 
positively.  Since  the  negative  resistance  is  low,  it 
is  essential  to  use  a  nearly  closed  structures  to 
reduce  the  loses  caused  by  radiation  of  energy. 
For  the  production  of  EM  oscillations,  the 
capacitance  of  the  diode  must  be  combined  with  a 
resonant  structure  having  the  proper  inductance  to 
resonate  at  the  desired  frequency,  and  having  a 
resistance  which  effectively  is  less  in  magnitude 
than  that  of  the  electron  beam.  Because  of  the 
low  losses  thus  required  of  the  device  operation 
the  properties  of  concentric  lines  and  of  tuned 


Fig.  3a  Ten-centimeter  diode  used  in  tests. 
Diode  ri  T2  x  h 

No.  24  1.270  0.635  0.203  1.870 
No.  37  1.220  0.635  0.105  1.870 
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cavities  offer  a  favorable  method  of 
approach.  The  experiment  was  carried 
out  using  Vo  =  300V  and  I  =  500mA. 
(no  gas  pressure  was  given).  Power 
output  and  efficiency  data  obtained 
during  the  experiment  are  shown  in  Fig. 
3b.  Power  output  of  a  few  tenths  of  a 
watt  at  efficiency  of  tenths  of  a  percent 
are  achieved.  Again,  there  is  no 
convincing  evidence  from  this 
experiment  that  identifies  the  physical 
mechanism  for  producing  the  EM 
radiation. 

3.  Microwave  Radiation  From 
Electron-Rich  Sheath-Plasma 

Resonance 

By  placing  a  small  size 
electrode  (with  a  positive  dc  bias), 
Fig.  3b.  Power  output  and  plate  efficiency  which  can  be  a  plane  disc,  sphere,  or 
for  diode  no.  24.  cylinder,  into  large  size  afterglow 

plasma,  microwave  radiation  with 
copious  harmonics  was  detected.  In  this  experiment7,  the  small  electrode  acts  as  an  rf 
antenna.  A  typical  discharge  can  be  divided  into  plasma  region  with  equal  electron  and 
ion  densities  and  an  electron-rich  plasma  sheath.  The  sheath-plasma  resonance  is  a  well- 
known  phenomenon  of  rf  probes  and  antennas  in  plasmas.  Basically,  an  electron-rich 
sheath  acts  like  a  capacitor  (diode)  while  a  field- free  cold  plasma  can  behave  inductively 
©  2 

(  s  =  1  -  <  0  for  ©  <  cop )  such  that  a  resonance  is  possible.  The  instability  is  driven 

by  a  negative  rf  sheath  resistance  associated  with  the  electron  inertia  effect  in  the  diode 
like  electron-rich  sheath.  The  maximum  values  which  the  negative  resistance  attains 
occur  in  the  neighborhood  of  transit  angles  given  by 

9  =  ©T  =  27m  +  ~,  n=l,2,3,....  (3) 

where  x  is  the  electron  transit  time  as  it  streams  through  the  sheath  region. 

Using  a  second  movable  floating  rf  probe,  the  field  amplitude  has  been  mapped 
versus  radical  distance  from  the  exciter  probe  (Fig.  4a).  It  shows  that  both  the 
fundamental  and  the  second  harmonic  exhibit  a  monotonic  amplitude  decay  in  the  near- 

c 

zone.  With  increasing  distance  (r  >  —  ~  1.6  cm)  the  fundamental  amplitude  truly 

(Dp 

vanishes  while  the  harmonic  remains  observable  throughout  the  plasma  and  beyond.  The 
result  clearly  indicates  a  highly  localized  field  pattern  of  the  nonradiating  fundamental 
mode  and  a  broader  pattern  of  the  propagating  second  harmonic.  This  also  confirms  that 
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the  plasma  behaves  inductively  (©  <  cop).  The 
observed  field  polarization  observed  (ref.  7)  is 
inconsistent  with  other  radiation  mechanism 
proposed  for  earlier  observations  of  EM  emissions 
from  positively  biased  wires  in  plasmas  such  as 
beam-plasma  instability6  and  orbitron 
mechanism5.  The  most  direct  confirmation  that 
the  above-mentioned  instability  mechanisms  are 
unable  to  explain  the  present  observation  is  shown 
in  Fig.  4b.  EM  harmonic  lines  are  observed  for 
all  basic  electrode  geometries:  plane,  spherical, 
and  cylindrical.  It  is  inconceivable  that  long- 
lived  electron  orbitron  or  counterstreaming 
electron  beams  arise  at  a  large  plane  electrode, 
whereas  sheath-plasma  instabilities  can  occur 
for  any  electrode  geometry.  The  observed 
instability  exhibits  threshold,  growth,  and 
saturation.  The  rf  voltage  and  power  saturate 
at  levels  small  compared  to  the  dc  voltage  and 
power,  respectively.  The  conversion 
efficiency  is  about  5x1  O'5.  This  low 
efficiency  may  be  attributed  to  that  the  small 
size  antenna  radiates  into  large  volume  plasma 
and  the  lack  of  resonator  to  confine  the  EM 
radiation. 


-JO  -5  0  5  » 

RADIAL  POSITION  Ar  (mm) 


III.  Numerical  Model  Fig.  4.  (a)  Radial  rf  amplitude  profile 

showing  a  broader  half- width  for  the 

The  efficiency  of  Terahertz  sources  propagating  harmonics  than  for  the 
based  on  the  DC  discharge  in  coaxial  evanescent  fundamental  emission, 

resonator  may  be  substantially  improved  from  (b)  Comparison  of  emission  spectra 
10'  provided  that  the  generation  mechanism  for  different  electrode  geometries 
can  be  identified.  In  order  to  fulfill  this  goal,  showing  the  existence  of  the  sheath- 

the  best  approach  is  to  build  a  particle  code  plasma  instability  for  all  shapes  and 

that  includes  electrostatic  and  electromagnetic  ruling  out  the  orbital  instability 
interactions  as  well  as  the  DC  discharge  mechanism, 

physics.  Since  the  electrostatic  field  plays  a 
very  important  role  in  determining  the 

electron  equilibrium  orbit,  we  will  start  with  building  an  electrostatic  particle  code 
bounded  by  coaxial  resonator.  The  B-spline  of  fourth  order  (degree  3)  is  used  for  solving 
the  two  dimensional  Poisson  equation  in  cylindrical  coordinates 


a2®  i  ao  a2® 

Sr2  +r  dr  +  8z2  ““4'tp  ’ 


(4) 
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with  the  boundary  conditions  of  ®(r,0)  =  <J>(r,L)  =  ®(rw,z)  =  0.  Expand  the  functions  O 
and  p  along  the  z-direction  in  terms  of  sine  functions 


<E>=  X  <£k(  r  )  sin  *j-  z  , 
k=l 


V  ,  \  ■  ** 

:  X  pk(r)sm-j^  z  , 

k=l 


where  Nz  is  the  number  of  grid  points  in  the  z-direction.  Substituting  (5)  into  (4)  yields 
1  1c7T  2 

Ok"  +  -  Ok  -  ( l"  )  Ok  =  -  47tpk  ,  (6) 

A  spline  can  be  expressed  in  the  following  form 


S(r)  =  XaiBi(r)  =  Ok(r)  ,  (7) 

i=i 

where  n  =  Nr+  2  (Nr  is  the  number  of  grid  points  in  the  radial  direction). 


The  detailed  property  of  B-spline  has  been  described  in  Ref.  [11].  In  a  coaxial  resonator, 
additional  boundary  conditions  at  the  center  conductor  O(r0,z)  =  Vo  must  be  imposed. 
With  this  applied  potential,  the  resulting  surface  charge  on  the  inner  conductor  is 
-Vo 

a= - - —  and  the  radial  electric  field  is 

47:  r0  ln(  f ) 

ro 

v  n 

Er(ro)  =  - - T-  =  -  Z  «i  B,'(ro) .  (8) 

r0  ln(  — )  i=i 

L0 


r/r0  r/r0 


Fig.  5  Comparison  of  (a)  potential  and  (b)  radial  electric  field  profiles  between 
theoretical  results  (dark  curves)  and  that  from  the  B-spline  calculations 
(dashed  red  curves). 
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For  the  case  of  infinite  coaxial  cylinders,  the  theoretical  profile  of 


Er(r)=y  and 


In(-) 
V(r)  =  Vo  T~ 

Hf) 

Lo 


(9) 


are  displayed  as  black  curves  in  Fig.  5.  To 
test  the  accuracy  of  the  B-spline  approach, 
results  from  the  B-spline  solution  are  shown 
as  red  curves  in  the  same  figure.  The 
potential  contour  obtained  from  the  B-spline 
solution  for  the  case  of  p  =  0  in  a  two 
dimensional  system  (r,z)  is  shown  in  Fig.  6. 

Since  the  maximum  energy  an 
electron  can  attain  is  in  the  neighborhood  of 
10  kV,  the  relativistic  equation  of  motion 
employed  in  our  particle  code  for  simulating 
microwave  devices  is  simplified  to  become 
non-relativistic.  In  order  to  allow  all 
possible  mechanisms  to  occur,  three 
dimensional  electron  orbital  motions  will  be 
allowed.  The  electron  equilibrium  orbit  can 
be  derived  from  the  Lagrangian  of  the 
energy  equation 


Fig.  6  Two  dimensional  (r,z)  potential 
contour  calculated  from  the  B- 
spline  solution. 


m  r  ,  dr  2  2  d0  2 , 

7  [(^)  +r  (*>  ] 


eV  (r) 


(10) 


The  cross-sectional  view  of  four  representative  electron  orbits  starts  with  different  initial 
condition  is  shown  in  Fig.  7.  Here  electrons  start  with  the  same  kinetic  energy  that  is 
equal  to  the  potential  energy  at  the  initial  position  but  at  different  angle  \p  with  respect  to 
the  radius  vector. 


The  electromagnetic  field  can  be  taken  into  account  by  using  the  same  technique 
as  we  have  employed  in  simulating  microwave  devices11.  The  scalar  function  appropriate 

to  the  lowest  E-mode  in  a  coaxial  guide  is  Coo  =  In  r  /  A  /  2tt  In  ~  .  Hence  the  field 

V  r0 

components  of  this  mode  are 

Er  =  Et(z) - 7= 


^lna 

r0 
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1 


B$  =  Bt(z) 


1 


(11) 


E(j,  =  Ez  =  Hr  =  Hz  =  0 

and  the  corresponding  Maxwell’s  equations  are 


5Er  dBj, 

irc(-  0z 


)  -  4tc  Jr 


(12) 


5Ba  SEr 

at  • 


where  Et(z)  and  Bt(z)  can  be 
expanded  by  using  Cosine  and  Sine 
functions  respectively  to  satisfy  their 
boundary  conditions.  This 

transverse  electromagnetic  mode  is 
the  dominant  or  principal  mode  in 
coaxial  guide.  If  the  higher  E-modes 
are  needed,  they  can  be  derived  from 

r 

the  scalar  functions  C,  =  Zm  (y,  ~ ) 

To  7 

m0  and  be  included  by  keeping 

olxl 

only  a  few  m  modes. 


Fig.  7  Example  of  electron  orbits  in  a  logarithmic 
potential  from  particle  simulations, 
v)/  =  (a)  20°,  (b)  40°,  (c)  60°,  (d)  80°. 
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